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ABSTRACT 

The analog computer has been used to obtain solutions to the diffusion 
equation for the case of a continuous point source aloft. The basic model 
used is that of a source located between an inversion and the ground vhere 
the -wind speed and eddy diffusivity are constant over the' diffusing region. 
Solutions were also obtained with several special features added to the basic 
model: an increase of lateral diffusivity over the vertical; representation 
of the eddy diffusivity as a function of plume dimensions or as a function of 
distance downstream; presence of wind direction shear in the diffusing layers; 
and the presence of large-scale eddies which cause the meandering of the 
plume . ' 

The solutions are presented as deviations from the reference solution^ 
that solution derived from the basic model^ to show graphically the effects 
of each parameter change. 
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1„ INTRODUCTION 

Turbulent diffusion from an elevated point source is of primary concern 
in air-pollution meteorology. Although such diffusion is the central compo- 
nent of a common problem^ it is by no means a well-understood phenomenon „ 
This results in part from the lack of a sufficiently general mathematical 
model with a complete set of solutions „ No matter how well diffusion is 
documented in the field ^^ it cannot be completely and compactly expressed 
and understood except in reference to such a model » But atmospheric diffu- 
sion is such an intricate mechanism that no one model is ever likely to 
satisfy all damandSo 

One model which is perhaps sufficiently general is the parabolic dif- 
fusion equation 
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and 



X = average concentration 

x^y^z = space coordinates 

u^v^w = average wind components 

Kj_ = coefficient of eddy diffusivity in the i-direction 

t = time 

1 



This model has tvo main shortcomings s it yields only average concentration, 
and a complete set of solutions has not been :fo.undo In fact, such analyti- 
cal solutions as are available apply only to special cases of limited appli- 

1 2 
cation o ^ Another deficiency is that the model implies an infinite diffus- 
ing velocity, but in cases where this is important, the hyperbolic diffusion 
equation can be used. In the case of steady- state diffusion, the velocity 
of diffusion is not important. 

The validity of the application of the parabolic diffusion equation to 
atmospheric diffusion has been questioned,-^ but the lack of adequate solutions 
has made verification difficult. The use of electronic computing machines 
makes solutions available for any type of parameter variation but here an- 
other difficulty arises o Since analytical solutions are not readily avail- 
able, the functional dependence of some of the parameters has not been de- 

h 
terminedo For example, Richardson has suggested that the coefficient of 

eddy diffusivity be a function of plume size, but the exact form of this 

function has not been indicated. 

The purpose of this paper is to show how the diffusion equation may be 
implemented on the electronic analog computer and to show the qualitative 
nature of the solutions. Means for obtaining improved solution accuracy are 
also indicated'o The reasons for using the analog computer in this applica- 
tion are that it is easy to program and that it presents the solution in 
graphical form. 

The problems chosen for solution:are based upon the case of diffusion 
in the steady- state from a point source located half-way between the surface 



of the ground and the base of an inversion. The special cases considered 
here were selected from the many problems of interest described else- 
where. " The methods outlined are perfectly general^ however^ and with an 
analog computer of sufficient capacity^ a similar analysis may be made for 
nearly any meteorological conditions Added to the basic model are several 
functions of the average wind components and the coefficients of eddy dif- 
fusivity so as to consider increased lateral over vertical diffusivity, the 
case of a plume passing from one turbulent regime to another^ the theoretical 
question of making the diffusivity a fun.ction of cloud size^ the case of 
wind direction shear with height in the diffusing region^ and the case of 
the meandering plume. The model is arranged so that the origin is at the 
source^ the x-axis oriented downwind^ the y-axis crosswind^ and the z-axis 
vertical The basic model includes no provision for variations in the wind 
speed with height or in the eddy diffusivity in the y- or z-directions. This 

Q 

can be readily provided'^ but the emphasis here is on special techniques. 

2, FINITE DIFFERENCE MODEL 

The model presented in Eq« (l) can be simplified by deleting the terms 
that do not apply^ to 

where the diffusion term in x was neglected in comparison to the transport 
term in x . The boundary conditions are 



as x^y^z -> 00 



K^ -^ -> as z -^ - -, + - f or X > 



where H is the height above ground of the base of the inversion layer « 
These conditions state that the concentration is infinite at the source 
and that the ground and the inversion are impervious to the diffusing sub- 
stance. The continuity condition is 




(5) 



for X > 

where Q is the rate of emission from the source. It is convenient to 
state (2) in nondimensional form: 
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where 

s = x/xo 

Y = y/H 
Z = z/H 

X = H^ u x/K^ 

V = Hv/Kz 

and the symbol X is reserved for later use. The constant Xq will be de- 
fined later. 



The mathematical model to be used^- Eq., {k) y must be put in a form suit- 
able for programming on the analog computer o Since the computer can integrate 
only -with respect to one independent variable., computer time, two of the in- 
dependent variables must be replaced with finite differences o The deriva- 
tives in Y : and ■ Z will be expressed as finite differences as follows: 

(5) 
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The expressions for the second derivatives in Y and Z are different be- 
cause, in all the cases studied here, uniform differences are taken in Z but 
nonuniform differences are taken in. Y o If uniform differences were taken 
both ways, the two expressions would be identical.. Equation {k) may now be 
stated as a set simultaneous, ordinary differential equations 2 
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where the subscripts m and n refer to discrete values in Y and Z , 
respectively. The physical model corresponding to Eq. (6) is shown in Fig. 
1. This is the particular model that will be used in the section on compu- 
tation in connection with the isotropic diffusion caseo It consists of an 
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Fig. 1. The 8U-cell^ finite difference model of isotropic diffusion 
from a point source. 



array of cells of rectangular cross section with width AY y height AZ ^ and 
with one end at x = and the other at x = 00,0 The cells are stacked 
vertically from the ground to the inversion base and horizontally from ■- ^ 
to + 00 o The number of cells in Y is finite since the end cells have width 
AY = 00 o The source cell is in the center of the array and is designated 0^0. 

All the properties of each cell are lumped in the center of the cross 
section so -that the model is an array of rods with one end at x = and the 
other at x = 00 , The source is no longer a point source but a vertical area 
source. In this models the continuity Eq, (j) holds for all values of x^ 
and by evaluating it at x = ^ we can define j^ . 



Q = \ \ H'^ u S Xo <aZ clY 



^1 ^^1 



The integration at x = is carried over the source cell only^ and in a 
given cell all the parameters in the integrand are constant*, Since the value 
of S at X = in the source set at unity^ 

t2 ~ 



Q = H^ u ^ (AY)^ (AZ)^ 



and therefore 



Q 



"" = u h2 (AY)^ (az)o ° 

The solutions that will be obtained with this model will be the average con- 
centration over the cross section of a cell as a function of its lengthy or 
distance downwind » 
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3. ANALOG MODELS 

In analog computation^ it is helpful to realize that. there exists a 
correspondence between the components of the circuit set up to solve Eq. 
(6) and the physical model shown in Fig, lo This can be shown by con- 
structing on paper a passive network analog as shown in Fig. 2. In this 
circuity there is a capacitor for each cell and a network of resistors con- 
necting them. In this analogy the voltage at each node is a function of 
time and corresponds to the concentration as a function of distance down- 
stream. The current into each capacitor is given by 

dV ^ ^ ^ dx 

1 = c ■ — = u Ay Az 

dt dx 

and c = u Ay Az 

y Ky Az ' ^ K^ Ay ' 

In this circuit the boundary conditions are implemented directly. 
Since there must be no flux through the ground or the inversion, the' counter- 
part of flux, current, is made zero by simply making no connections across 
the boundaries. The initial concentration (voltage) in the source cell is 
provided by disconnecting the capacitor for that cell, charging it to ar- 
bitrary voltage, and switching it back into the network at x (time) 
equals zero. Since the end cells in y are infinite, the concentration 
(voltage) is always zero there so these nodes are connected to ground 

Provision for variation of the parameters K , u , Ay , and Az as 
a function of y and z is made by assigning appropriate values to the re- 
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Fig. 2. The passive network analog of diffusion from a point source. 
The center cell is the source cell and its capacitor is initially 
changed to represent the source concentration. 



sistors and capacitors involvedo If any of these parameters are a function 
of X ^ then they must be made to change in time« This is not readily ac- 
complished and constitutes one of the disadvantages of the passive network. 
In general^ the passive network is much less flexible than the active cir- 
cuit. 

An electronic analog circuit will be shown later for each case. They 
are constructed by assigning to each component a mathematical operation and 
these are connected so as to satisfy the governing set of differential equa- 
tions. Identical circuits can be derived by simulating the passive network 
analog; thus each circuit has a dual nature. It is a series of mathematical 
operations designed to solve an equation and it is also an analogy of the 
system being studied where components and voltages represent physical parts 
of the atmospheric diffusion model. Both the active and the passive analogs 
are physical systems designed to behave like another^ idealized physical sys- 
tem. The electronic analog circuits will be explored as they appear. 

k. COMPUTATION AND RESULTS 

Since each of the cases considered is quite different^ the form of Eq. 
(6) appropriate to each will be given along with the computing circuit and 
the results. 

A. ISOTROPIC DIFFUSION WITH AN 8^-CELL MODEL 

The first case is diffusion from the point source in isotropic condi- 
tions^ that is^ the wind speed and eddy diffusivity are constant over the 

region. The 8^-cell model used is shown in Fig. 1„ It is symmetrical about 
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the coordinate axis so only the cells in the first quadrant need be simu- 
lated and of these, the ones in the outer column, m = 3 , are always zero 
so that there are only ik active cells in the simulation » In this case, Eq, 
(6) becomes 



dS 




1 


^+l,n ■ ^,n 


^,n " ^-l,n 


dX 


m,n 


196 (AY)jj^ 


2 


2 



Sm,vtM-l ~ 2Sni n + Sm n-1 
+ — ' ^ — 7^ ^ . 



196 (AZ)' 



(7) 



where X = 196 x and the boundary conditions are 



^l,n 


= °-l,n 


^,1 


= ^,-1 


Si, 1(0) 


= 1.00 


Sj^n 


= 


Sm,7 


= %8 



The factor 196 was inserted to make the coefficients more convenient o The 
cell boundaries given in Table I show that the spacing in the vertical is 
linear whereas the spacing in the horizontal is exponential and the source 
cell has a square cross section. 
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TABLE I 

CEIL BOUNDARIES IN THE HORIZONTAL AND THE VERTICAL 
FOR THE FIRST QUADRANT OF THE 8i(-CELL MODEL 



Y = y/H 


Z = z/H 








0.07li+ 


O.OTlil- 


0.195^ 


0,1^28 


00 


0.2142 




0.2856 




0.3570 




0.i^281+ 




0,5000 



Since AZ is a constant for each cell, Eq. (7) can be conveniently ex- 
pressed as two sets of difference-differential equations for the colimins 
used. 



dS 
dX 

dS 
dX 



^ = 0.8228(S2^n - Si,n) + (Si,n+1 " 2Si^n + S^^n-l) 
l,n 

= 0.1.937 Si^„ - 0.7^10 Sg^^ 

The computer circuit for Eq. (8) which is shown in Fig. 5 has one ampli- 
fier for each of the 1^ differential equations represented by (8). The top 
row of amplifiers simulates the first coliJinn of cells and the top equation in 
(8)o The bottom equation is solved by the lower row of amplifiers. The Sj_ j^ 
amplifier in the top row represents the source cell and has an initial condition 
of S = 1.00 at computer time (distance in X) t = . The basic unit of 
computation is the integrating amplifier as shown in Fig. k. for a typical cell. 
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Fig. 4. Typical cell from the circuit of Fig. 5 sho-wing the connections 
to the integrating amplifier. . Nonunity equation coefficients are set on 
the coefficient potentiometers provided. - 
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Solution curves for the first rov of cells are shown in Figs, 5 and 6. 

Bo ISOTROPIC DIFFUSION WITH A 56-CELL MODEL 

To facilitate the simulation^ the number of cells can be reduced. In 
this case^ only k cells in Y are to be used instead of six. The purpose 
of this reduction is that more involved problems to be treated are more readily 
handled with a simpler basic model « Not so much equipment is needed in this 
case but the reduction in the number of cells imposes a severe accuracy pen- 
alty. Before introducing the more interesting cases ^ the magnitude of the 
error introduced can be shown by solving the same problem employing the 
same techniques. 

Cell boundaries in Y for the 56-cell model are Y = , O.O71U ^ 00 ; 
while the cell boundaries in Z are as stated in Table I. In the first 
quadrant^ there are still 7 active cells in Z but only 1 in Y ; thus there 
is 1 active column and Eqo (8) simplifies to 

- 0.7220 Si^n + (Si^n+1 - 2Si^n + Si^n-l) (9) 



dS 
dX 



where X = 196 x and the houndary conditions are 
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The computer circuit for Eq. (9) is shown in Fig. 7° As in Fig. 3^ there is 
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Fig. 5- Concentration plots for each cell as a function of distance down- 
wind from a point source using the 8^4— cell model in isotropic conditions. 
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Fig. 6. Concentration plots for each cell as a function of distance down- 
wind from a point source using the 81+-cell model in isotropic conditions. 
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one integrating amplifier for each active cello The solutions for the 
source cell obtained from the 56-cell and the 84-cell models are compared 
in Figo 8. Evidently the approximations are not bad for small values of 
X o For all the following work, the absolute error for values of X > 1 is 
large, the objective is to shov deviations from the reference solutions in 
each case. 

C o NONISOTROPIC DIFFUSION 

Changing the ratio of the eddy diffusivity in the Y-direction to that 
in the Z-direction produces a very marked alteration in the solutions o As 
an example, set Ky = 10 K^^ in the 56-cell model o Equation (9) may now be 
written 



dX 



^ ^ = - 0o7220 Si^^ + — ^— ^ — ^— (10) 



where X = 196O x and the boundary conditions are 



l,n -"l^n 

S = S 

-i- • JL JL 4 ^ J« 

S (0) = loOO 

So =0 
2,n 

^1,7 " ^1,8 



The computer circuit for Eqo (lO) is identical to the one shown in Fig, 7 

except that the coefficient potentiometers are set to the values indicated 

in Eq„ (lo)o Solutions for the cases Ky = K^ and Ky = 10 K are plotted 

in Figo 9» It shows that increasing the lateral diffusivity in this model 
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Fig. 8. Comparison of concentration plots for the source cell using the 
8k- and the 56-cell models with isotropic diffusion conditions. 
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Fig. 9. Comparison of concentration plots for some cells using the 56-cell 
model. This plot shows the changes produced in the solution by increasing 
the lateral diffusivity by 10 times over the vertical diffusivity. 
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has a pronounced effect. The magnitude of the effect is enhanced by the 
limited vertical diffusion due to the presence of the inversion. 

D. EDDY DIFFUSIVITY RELATED TO PLUME DIMENSIONS 

Plume dimensions may influence the magnitude of the effective eddy dif- 
fusivity acting on the plume since^ as the cloud grows ^ larger eddies are 
able to operate to diffuse the materials. The eddy diffusivity might be 
made a function of the plume size and of the power spectrum of eddies pres- 
ent. In the usual parabolic diffusion equation^ the plume has no easily de- 
fined boundaries except those imposed upon it^ such as the surface of the 
ground. Therefore the diffusivity cannot be related directly to plume size 
but must be specified in some arbitrary manner. For example^ it could be 
made proportional to the standard deviation of the concentration distribution 
taken perpendicular to the plume axis. If this distribution is normal, the 
standard deviation is inversely proportional to the central maxima or axial 
concentration, so that the diffusivity can be set inversely proportional to 
the axial plume concentration. Since the power spectra would never be con- 
tinuous indefinitely, the diffusivity could not become infinite as the axial 
concentration went to zero but would reach some finite value and remain con- 
stant. In the finite difference model, the source is finite so that the dif- 
fusivity wolild increase from one small finite value to a larger one inversely 
as the concentration in the cell containing the source. Represent the dif- 
fusivity as 0(Sq_)K where K is the diffusivity for infinite plume dimen- 
sions and ^{Sj) is a function of the axial concentration such that 

< 0(Si) < 1.00 . 
22 



For example, set 

0(S-l) = 1 - O.T S^^-L (11) 

since 1.00 > S^, i ^ o 

Setting Ky = Kg , Eq. (9) "becomes 



dS 
dX 



= i25(Si)(Si^^^.l - 23^^^^ + Si^^.i) - 0.7220 0(Si) S^^^ (l2) 

x^ n 



with the same boundary conditions as before « The computing circuity shown 
partially in Fig. 10^ is similar to that of Pig*. 7 ^^ principle. This one uses 
5 gainplifiers^ one integrator and two summers^ per cell to isolate all the in- 
puts to the integrating amplifier. With this configuration j, it is only neces- 
sary to introduce the modifying parameter 0(S]_) once per cell. 

Since 0(S]_) is a function of the axial concentration, its form is de- 
termined by a feedback look defined by Eqs„ (ll) and (13). The function ^(Sj^) 
is generated as the solution proceeds.. 

The computing circuits used in this problem are virtually identical with 
those of the next section. Therefore the description of the additional com- 
ponents i*equired for the generation of 0(S-^) and the presentation of the 
solutions will be deferred to that section o 

E. DIFFUSION IN TRANSITIONAL STATES 

A problem of special interest arises when the trajectory of the diffusing 
plijime crosses the boundary from one turbulent regime to another. The physical 
model consists of a source located over land not far from the shore of a lake. 
The wind blows from the land out over the lake and the mechanical and perhaps 
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the thermal turbulence decrease markedly as it does so. The associated dif- 
fusion is referred to in this research as diffusion in a transitional state; 
the inversion breakup fumigation is the result of a corresponding time transi- 
tion in diffusion. 

In the mathematical models a change in the eddy diffusivity accompanies 
the transition from one turbulent regime to the other. This vill be implemented 
by applying another modifying parameter 0(X) to the constant K , Then 0(X) 
will be unity in the first regime and fall to some lesser value in the second. 
For the present preliminary study^ this function was simulated directly on the 
computer to conform to an intuitive impression of what it should be without 
using a mathematical expression. It is possible to contrive a mathematical 
description^ which is: 

0(x) = 0,52 + 0A8([]l - h(X - Xi)3 

+ cos(X - Xi) \ji{X - Xi) - h(X - Xi - It)]} (15) 

where h(a) is the unit step function which has the following properties: 

h(a) = ^ a < 

= 1 , a > , 

and X2^ is the value of X where the turbulence starts to decrease „ 

The diffusion eq.uation is just Eq. (l2) with 0(Sj_) replaced with 0(X) . 

The computing circuit is also that of Fig, 10 and the circuit for generation 

of 0(S2_) and 0(X) is shown in Figo 11. 

The plysical problems discussed in this section and in the previous one 
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MECHANICAL LINKAGE TO 
FOLLOW-UP POTENTIOMETERS 

Fig. 10. Computing circuit for the first quadrant of the 56-cell model which 
permits the inclusion of special forms of the coefficient of eddy diffusivity. 
This configuration permitted the generation of solutions where the diffusion 
coefficients were a function of plume size and in the case of diffusion in 
transitional states. The dotted lines indicate mechanical connections. 
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are different despite a some^^hat deceptive similarity in the governing 
equation used. The differences between the problems are indicated in the 
form of the modifying parameter » It is only because of the similar way 
in which is implemented that the problems are handled together » 

It is possible to introduce both effects in the problems o In this case^ 
a new parameter will be defined such that 0(S]_^X) = 0(Si) 0(X) „ Figure 11 
shows how each of three functions is generated o When switch 1 is on -100 
and switch 2 is closed^ the output of amplifier k is 0(S-]_) as defined by 
Eq. (ll) and it passes through servomultiplier 1 unchanged to control servo- 
multiplier 2. When switch 1 is on 0(X) and switch 2 is open, the function 
0(X) , generated by amplifiers 1^ 2, and 5^ controls servomultiplier 2, The 
function 0(S]_,X) is generated when switch 1 is on 0(X) and switch 2 is 
closed. The fimctions ^{Sj) ^ 0(X) and 0(Sx,X) are plotted in Fig. 12 as 
they were generated by the computer o 

Solutions for S-j_ q^ , ^12 ^^^ ^1 5 ^^^ presented in FigSo 15^ 1^, and 
15^ respectively, for k cases 2 = 1, = 0(S-l) , = 0(X) and = 0(S3_,X) 
The case 0=1 is the reference solution „ Each modifying parameter acts to 
change the scale in X , that is, the new solutions can be found by transla- 
tion in X from the reference solution. This is shown best in Figs„ 1-^ and 
15 where the maxima have been translated and expanded but the value of the 
maxima remains the same. The first modifying parameter 0(S]_) is signifi- 
cantly different from unity only for small values of X so the greatest :■: 
effect on the solutions should occur there o Each solution shows large trans- 
lation for small values of X » The net effect of coupling the eddy diffu- 

sivity to cloud dimensions is to increase the concentration values for small X 
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Fig, 12. Plot of the modifying parameters of the coefficient of eddy dif- 
fusivity. 0(S]_) is a function of plime size while 0(X) represents trans- 
ition from one turbulent regime to another. 0(S2_,X) is a simple multipli- 
cative combination of 0(X ) and 0(X) . 
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2 3 4 5 6 7 

DISTANCE DOWNWIND X 

Fig. 13. Plot of the concentrations in the source cell as a function of 
the " modifying parameter used. 
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Fig. lU. Plot of the concentration in cell 1, 2 as a function of the 
modifying parameter used. 
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2 3 4 5 

DISTANCE DOWNWIND X 



Fig. 15. Plot of the concentration in cell 1, 5 as a function of the 
modifying parameter used. 



31 



along the pliame axis and to delay the occurrence of the maxima off the axis. 
The solutions for 0(l) and 0(Sj_) do not converge for large X but reach 
a constant translation in X . But since the gradient is small for large X ^ 
the absolute value of the differences between the solutions is smalls 

Since the function 0(X) can be determined independently^ the solution 
could be found from the reference solution by plotting it on a nev scale de- 
fined by 



= r 0(x) 



•X 
X-^ = 1 0(X) dX 




The observed effect corresponds to this point of view in that 0(X) has no 
effect on the solution for X < 1 but acts to translate the solution for 
X >1. 

F. WIND DIRECTION SHEAR 

The presence of wind direction shear in the diffusing layers causes 
the plume cross section to be skewed, and the amount of skewness can be de- 
termined by inserting terms in the diffusion equation for cross-wind components « 
Equation (6) contains the nondimensional parameter V = H v/K for these 
cross-wind components. Since the diffusion pattern is not symmetrical, it 
will not suffice to simulate only the first quadrant. In the interest of 
simplicity, a new 25-cell model is used as shown in Figo l6o The source 
cell is 0,5 and there are 15 active cells. Figure 17 shows V(Z) , which is 
an exponential function of height, in the finite difference model in relation 
to the downwind component H u/K . The exponential function was chosen as an 
approximation to the lower part of the Ekman spiral. Equation (6) may be 
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written 



dS 
dX 



1 



25(AY)in 



(AY)m+l ■ (AY)m.l 
2 2 



/^Y) L^+l^"^ " ^-l>nj "*" L^m^n+l - 2Sm^n + Sjji^n.iJ 



50 (AY) 



(li^) 



where X = 25 x and the boundary conditions are 



S.2,n 


= 


S2,n 


= 


^,0 


= ^,1 


Sm,5 


= ^,6 


So,5(0) 


= 1.00 



The computing curcuit is shown in Fig, l8 and the solutions are in Figs. 
20^ 21^ and 22 « Each solution is compared to the reference solution, where 
V(Z) = . We can compare the concentrations in the cells for two cases at 
some distance downstream, say at X = 1 ."^ The results, as shown in Fig. 19, 
indicate that the wind direction shear produces a marked skewness in the 
plume. The skewed plume is to be expected but a more interesting question 
is whether the average concentrations in the plume are reduced because of the 
skewness. Figures 20, 21, and 22 indicate that the concentrations are reduced , 



■'^In making this comparison, note that changes in wind direction were provided 
by adding v components while holding u constant which increased the mag- 
nite of the total vector by an average 1.6^ o Since X is a function of u , 
we can compensate for this by evaluating the concentrations when V ;^ at 
X = 1 and when V = at X = l„0l6 . 
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(a) WITHOUT WIND SHEAR 




(b) WITH WIND DIRECTION SHEAR 
X = l.00 

Fig. 19. Cross section of the plume at X = 1 showing lines of constant 
concentration in arbitrary units. This figure shows the reference solution 
(a) and solution in the presence of wind direction shear with height (b) . 
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Fig. 20. Comparison of the concentration plots for cells to the right and 
the left of the source cell in the presence of wind shear. The solid line 
is the reference solution for both cells. 
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Fig. 21. Concentration plots for some cells in the 25-cell model in the 
presence of wind shear. The solid lines are reference solutions. 
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Fig. 22. Concentration plots for some cells in the 25- cell model in the 
presence of wind shear. The solid lines are reference solutions. 
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by about h$ in the center of the plijme at X = 1 • ■ The Integrated concentra- 
tion over the plume at X = 1 was reduced about 2^. Thus it seems valid to 
conclude that wind direction shear does increase the diffusing properties of 
the atmosphere. It is planned to analyze the influence of wind direction shear 
in concentrations downwind from a continuous point source using an analog com- 
puter of much greater capacity. 

G. LARGE-SCALE EDDIES 

The presence of large-scale eddies or^ since we are considering the 
steady state, a standing wave in the horizontal wind components can be repre- 
sented as a sinusoid v = a sin bx . A stationary plume might be an example 
of this kind of process. The small-scale .eddies . which operate to diffuse 
the plume are represented by the coefficient of eddy diffusivity, whereas 
the large-scale eddies act only to transport the plume in the y-direction. 

Equipment limitations at the time of the research dictated the choice 
of a more coarse finite difference model as shown in Fig. 25 « The mathem- 
atical model, based on Eq. (6), is 



dS 
dX 



'- ^m+l,n + ^m-l,n "*" '^m,n+l "^ ^,n-l " ^^m^n 
m,n 



" 18(AY) C^^"*-!^^ " Sm-i^n] (15) 

/ N Ha 
where V(X) = "K" sin BX and the boundary conditions are 



^2,n - ^-2,n = ^ 
So^5(0) = 1.00 
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and AY = AZ for the active cells. If we set 5A = Ha(AZ)/2K , (15) becomes 



dS I 

^ lm,n 10 



= Y5 [^^.^-1 " ^^^^^ + (1 - 5A sin BX) S^^.1^^ 



+ (1 + 5A sin BX) Sin-i^n] (16) 

where X = 90 x . The computing circuit is shown in Figs. 2k and 25, and a 
typical cell is shown in Fig. 26, The solutions are plotted in Figs. 27^ 28, 
and 29 along with the reference solution. 

This system exhibits a phenomena like resonance in a forced second^ order 
system. It is not pure resonance but the sytem does respond markedly to 
changes in frequency of the cross-wind component. Changes in amplitude do 
not affect the qualitative nature of the solution. 

The plume meander lags behind the cross-wind components. The point at 
which the plume axis crosses the middle cell can be found by noting the point 
where the concentration, of opposite pairs of cells is equal. This occurs for 
X = 1.506, 2»khQy and 5-551 , while the cross-wind goes to zero at X = 1,0^1-7, 
2,09^^ and 5-1^1- The phase lage is not constant for the above three points 
but is 78 • 9 deg, 6O.9 deg, and 67 ol deg, respectively. 

It is of interest to compare the concentrations in the plume, for samp- 
ling periods short enough to avoid plume shifting, when it meanders with when 
it does not. At an arbitrary point downwind, at X = 1.6 , the total concen- 
traion reported in the case with no meander was l6^ higher than in the case 
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Fig. 27. Plot of the concentrations for cells 0,2, 0,5, and 0,4 in the 
presence of large-scale eddies. The solid line is the reference solution. 
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Fig. 28. Plot of the concentration for cells -1,5 and +1,5 in the presence 
of large-scale eddies. The solid line is the reference solution. 
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Fig. 29. Plot of the concentrations for cells -1,2; +1,2; -1.1+; and +1 h 
in the presence of large-scale eddies. The solid line is the reference' 
solution. 
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vhere the plume meanders. As shown in the appendix^ this decrease in concen- 
tration in the meandering plume is due entirely to the longer path length 
followed by the pl\:ime. 

If one knew the total concentration across a straight plume as a function 
of distance downwind one could then compute the long-term concentration at a 
station subjected to a plume in meandering but otherwise identical conditions » 
The calculations may be so tedious as to be impractical since they would in- 
volve computing the length of the trajectories over the given interval^ ad- 
justing the concentration and integrating them. 

It would still be of interest to perform, this calculation on the computer 
using random large-scale eddies^ and it is intended that this topic will be 
the subject of a future report. 

H. A POSSIBLE COMBINATION OF TWO OF THE MODELS 

Two of the models presented above are related to some extent » The model 
in which eddy diffusivity increases with the plume dimensions can be related 
to that including large-scale eddies » In the first case^ the complete spec- 
trum of eddies can be included and the effect on the plume partially deter- 
mined. The solution is partially complete in that^ while the diffusing ef- 
fects of all the eddies are included^ the transport due to eddies much larger 
than the plijme dimensions is neglected. The model simply ignores eddies of 
a given scale until the pl\JU3ie dimensions grow to that scale « The transport 
phenomenon is shown in the latter model where the large-scale eddies are 
treated as a wind component, but here they contribute nothing to the diffusion 
process no matter how large the cloud becomes. 
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A combination of the two models might be an improvement. Let the wind 
components represent eddies larger than the pliome while all smaller eddies 
act through the eddy diffusivity. The wind components would become smaller 
with distance downstream whereas the eddy diffusivity would increase. The 
wind components incorporated in the model are defined as averages which 
raises the question of how much they should be allowed to fluctuate or to 
possess time - or distance-variant properties. An averaged wind speed could 
be defined as the average or mean value obtained by the use of a sampling time 
just large enough to exclude all the turbulent eddies which contribute to 
the diffusion process • If we allow the assiomption that the scale of eddies 
involved depends upon the scale of the plume ^ then the sampling time in the 
above definition would depend upon the plume size also. The wind- speed 
components could be represented as a Fourier series whose higher-frequency 
components are gradually eliminated as the plume grows ^ and conversely^ the 
coefficient of eddy diffusivity would increase as more turbulent eddies were 
added to its domain. This is a rather involved model even for analog computa- 
tion and no attempt, as yet has been made to set it up. 

5- SUMMARY AND CONCLUSIONS 

As the model has become more complicated, it has been necessary to use 
fewer cells and to tolerate greater error , The amount of equipment needed 
is proportional to the number of cells and to the complexity of the model, 
and all the research described here was done on one machine with a very 
limited amount of computing equipment The direct way to reduce error, then, 
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is to use a larger machine "which "would permit simulation of more cells. The 
errors are proportional to the fourth partial space derivative so that the 
error can be reduced by solving for a function related to the solution which 
possesses a smaller fourth derivative. Such a function could be found by per- 
turbation techniques. 

The parabolic diffusion equation has been shown to be of general utility 
in the variety of problems covered. It is a model which is well adapted to 
analog computation. With sufficient equipment the actual programming seldom 
presents any difficulty even in the most involved problems discussed here. 

The trapping situation was used throughout this report to limit the 
scope of the research^ and to facilitate handling the problems on a computer 
of given size. Extension of these methods to nontrapping situations can be 
readily done using the technique reported earlier. 
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APPENDIX 

In the case where the meandering plume was treated, there were two pos- 
sible mechanisms which could affect the diffusion of the plume. The one noted 
in the main body of the report was the increased path length of the trajectory. 
The other is that large-scale eddies might change the shape of the plume, 
thereby influencing the diffusion process « An attempt was made to evaluate 
this mechanism by correcting for the first one. 

When there is a cross-wind which causes ^ meandering, the magnitude of the 
total wind vector and the path length along the plume axis, as shown in Fig. 
50, has been increased. To compensate for this, the concentration for no 
cross-wind should be taken at a point further downstream in inverse propor- 
tion to the relative decrease in the magnitude of the wind vector. Taking the 
first zero crossing, at X^^ = I.506 , the compensating factor is 



X^ = X; 



2^ VI 4- sin^ BX 



where 

'X 
sin"^ BX = — ■ 



1 r^ 2 1 r n 

— \ sin^ BX dX = -^^ [BX^ - sin BX^ cox BXjJ = 0. 1+790 

1 Jq 1 

Then X^ = 1.506/17^79 = 1,851 » The following table gives the concentration 
in both cases. 
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TABLE A-I 

COMPARISON OF PLUME CONCENTRATIONS 
FOR THE MEANDERING PLUME AND THE STRAIGHT PLUME 



Cell 


Straight Plume 

X = 1.851 


Meandering Plume 

X = 1.506 


0,2 

0,3 

0,4 


0.0595 
0.0596 
0.0595 


0.0599 

.0617 

0.0599 


-1,2 

-1,3 
-1,4 


0.0412 
0.0422 
0.0412 


o.o4o6 

0.0420 
0.o4o6 


+1,2 

+1,3 
+1,4 


0.0412 
0.0422 
0.0412 


0.o4o6 
0.0420 
0.o4o6 


Total concentration 


0.4278 


0.4279 



Evidently the plume shape is a little different since the individual cells do 
not report identical concentrations^ but the sume of the concentrations for all 
the cells is not significantly different. 
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